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Abstract 

The  primary  objective  of  the  study  is  to  reproduce  the  dynamic  formation,  development  and 
burst  of  the  vortex  that  forms  on  the  leeward  side  of  a  70°-swept  delta  wing  pitching  in  sinusoidal 
oscillations  around  its  40%  chord  point.  Emphasis  is  put  on  the  prediction  of  the  hysteresis  loops 
that  form  in  the  aerodynamic  loads  history  but  the  visualization  of  the  unsteady  flow  around  the 
wing  is  also  introduced.  The  freestream  Mach  number  for  all  the  computations  (both  Euler  and 
Navier-Stokes)  is  set  to  0.2  while  the  mean  angle  of  attack  is  changed  between  12°,  22°  and  38°  for 
oscillations  with  a  semi-amplitude  varying  between  3°  and  18°.  The  moving  grid  approach 
implemented  in  the  flow  solver  NSMB  is  used  to  reproduce  the  motion  of  the  wing  and  the  time- 
dependent  flow  equations  are  integrated  with  the  dual  timestepping  method.  In  addition,  a  method 
for  co-processing  the  simulation  and  the  visualization  of  the  flowfield  is  presented  using  flow 
feature  extraction  to  reduce  the  amount  of  data  saved  in  this  large  data  storage  requiring  test  case. 

Nomenclature 


c 

Root  chord  length 

Cd?  Cn,  Cm 

Drag,  normal  force  and  pitching  moment  coefficients 

f 

Oscillation  frequency  in  Hz 

k 

Reduced  frequency  coc/2Uoo 

Moo 

Freestream  Mach  number 

T 

Period  of  one  oscillation  in  s 

t 

Time 

Uoo 

Freestream  velocity 

x,  y,  z 

Rectangular  wing-fixed  coordinates  (chordwise,  spanwise  and  normal  directions 
respectively) 

a 

Instantaneous  angle  of  attack 

ocm 

Mean  angle  of  attack 

Aa 

Semi-amplitude  of  the  oscillation 

CO 

Oscillation  frequency  in  rad/s 

Introduction 

Current  and  future  fighter  airplanes  are  designed  to  transiently  fly  at  very  high  angles  of  attack. 
This  super-manoeuvrability  is  now  required  for  dog  fight  manoeuvres  like  the  "Cobra  manoeuvre" 
for  instance.  These  new  types  of  manoeuvres  involve  high  pitch  rates  and  flight  at  incidences 
beyond  the  static  stall  angle  of  attack.  The  aerodynamics  around  a  manoeuvring  fighter  is  very 
complex  but  needs  to  be  understood  in  order  to  optimize  the  airplane  design.  The  good 
aerodynamic  performances  of  the  fighters  in  high  angle  of  attack  manoeuvres  are  often  reached 
thanks  to  the  use  of  a  delta  wing.  The  aerodynamics  of  these  wings  at  incidence  is  characterized  by 
the  formation  of  a  pair  of  strong  leading-edge  vortices  on  the  leeward  side  of  the  wing.  The 
vortices  create  a  zone  of  high  flow  velocity  with  low  surface  pressures  on  the  wing  producing  an 
additional  lift  compared  to  classical  rectangular  shapes.  At  very  high  incidence,  however,  the 
vortex  flow  breaks  down,  a  zone  of  recirculation  with  a  turbulent  wake  appears  and  the  lift 
decreases.  During  manoeuvres,  the  same  phenomena  occur  but  the  flow  has  to  adapt  to  the  moving 
planform  and  thus  time-lags  are  observed  in  the  dynamic  response.  The  prediction  and  the 
understanding  of  the  hysteresis  loops  that  form  a  history  of  the  aerodynamic  forces,  for  example, 
can  be  very  useful  for  the  designer  to  extend  the  limit  of  the  flight  envelope  during  manoeuvres. 


Paper  presented  at  the  RTO  AVT  Symposium  on  ‘ Advanced  Flow  Management:  Part  A  -  Vortex  Flows  and 
High  Angle  of  Attack  for  Military  Vehicles ”,  held  in  Loen,  Norway,  7-11  May  2001,  and  published  in  RTO -MP -069(1). 
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That  is  why  at  the  beginning  of  the  80’s,  the  interest  in  aerodynamic  studies  of  moving  delta 
wings  became  important.  This  led  to  numerous  wind  tunnel  experiments  about  pitching  or  rolling 
wings  and  numerical  studies  of  the  phenomenon  appeared  in  the  90's.  The  scope  of  the  present 
study  is  limited  to  pitching  delta  wings  and  rolling  wings  are  not  assessed.  One  of  the  first 
experiments  on  pitching  delta  wings  is  the  towing  tank  experiment  made  by  Gad-el-Hak  and  Ho7. 
The  paper  has  become  a  classical  reference  as  it  sets  the  bases  for  the  following  studies  on 
sinusoidally  pitched  delta  wings.  The  first  part  of  Ashley  et  al.'s  survey2  gives  a  more  complete  list 
of  the  other  experiments  and  theoretical  studies  made  in  the  late  80's. 

Concerning  70°-swept  delta  wings  which  is  the  geometry  chosen  for  the  present  study,  quite  a 
few  references  are  available.  Soltani  and  Bragg18’19’20,  Thompson  et  al.21’22,  LeMay  et  al.16  as  well 
as  Kowal  and  Vakili15  looked  at  the  breakdown  position,  loads  history  and  pressure  measurements 
for  pitching  wings  in  wind  tunnel  experiments.  Brandon3  reports  experiments  conducted  at  the 
NASA  Langley  Research  Center  on  pitched  flat  plate  wings  with  sweep  angles  0°,  45°  and  70°  as 
well  as  on  a  F- 18  configuration.  It  is  important  to  mention  here  that  even  if  the  sweep  angle  is  the 
same,  the  geometries  of  the  models  presented  in  these  papers  are  not  exactly  identical  (single  or 
double  bevel,  presence  of  a  body  to  fix  the  mounting  system,  ...)  which  can  have  a  major  impact 
on  the  results,  as  shown  by  Ericsson  and  Beyers4.  The  paper  by  Hummel  and  Loeser11,  although 
for  a  different  geometry  (65°  cropped  delta  wing),  has  been  a  good  source  of  inspiration  for  the 
present  paper,  especially  the  aerodynamic  forces  results. 

On  the  numerical  side,  few  studies  have  pitching  delta  wings  for  subject  and  even  fewer  deal 
with  a  70°-swept  wing.  The  65°  cropped  delta  wing  experimentally  studied  by  Hummel  and 
Loeser11,  is  also  the  geometry  chosen  by  Fritz5’6  and  Arthur  et  al.1  for  numerical  simulations.  A  75° 
delta  wing  in  a  ramp-type  pitch  manoeuvre  is  simulated  by  Visbal  and  Gordnier23,24  solving  the 
Navier-Stokes  equations.  Visbal25  also  described  the  three-dimensional  instantaneous  structure  of 
the  flowfield  using  the  critical-point  theory  for  the  same  ramp-pitched  wing.  Kandil  and  Chuang12 
reports  the  surface  pressure  and  static  pressure  contours  given  by  an  Euler  simulation  for  a  75° 
delta  wing  pitched  around  20.5°  with  an  amplitude  of  2°.  Kandil  and  Abdelhamid14  later  performed 
Navier-Stokes  simulations  around  a  76°-swept  delta  wing  pitched  in  a  ramp  motion  up  to  90°.  A 
transonic  test  case  is  also  computed  by  Kandil  and  Kandil13  for  a  65°  delta  wing  pitched  around  20° 
with  4°-amplitude  oscillations.  The  complex  shock-vortex  interaction  is  shown  by  static  pressure 
contours  on  the  wing  and  plane  of  symmetry  as  well  as  velocity  vectors. 

The  present  paper  reports  the  state  of  the  work  started  at  the  Department  of  Aeronautics  of  the 
Royal  Institute  of  Technology  on  the  simulation  of  pitching  delta  wings.  A  previous  paper  by  Le 
Moigne17  gives  more  details  on  the  numerical  methods  used  in  the  study  and  presents  the  first 
results  obtained  for  (mainly)  Euler  calculations.  The  main  results  from  this  previous  paper  are 
recalled  here  and  new  turbulent  calculations  are  presented.  New  visualizations  of  the  unsteady 
flowfield  are  also  included  in  this  paper.  Early  in  the  beginning  of  these  simulations  of  a  pitching 
delta  wing,  arose  the  concern  of  the  size  of  the  result  files  saved  to  disk.  It  is  indeed  both  a  3D  test 
case  that  requires  fine  meshes  with  around  one  million  cells  and  an  unsteady  case  which 
necessitates  the  storage  of  several  time  accurate  solutions  for  its  analysis.  That  is  why  the  idea, 
developed  by  Haimes  and  Jordan10,  of  having  both  the  flow  solver  and  the  visualization  tool 
working  together  in  order  to  save  to  disk  only  the  relevant  part  of  the  flowfield,  has  been 
examined.  The  flow  solver  NSMB  used  in  the  study  has  been  slightly  modified  to  communicate 
with  the  visualization  software  EnSight  so  as  to  let  the  user  choose  what  part  of  the  flowfield  to 
extract  from  the  solution  and  save  permanently  to  disk.  Using  this  method  the  size  of  the  solution 
files  can  be  reduced  by  several  orders  of  magnitude. 

In  the  paper,  the  main  characteristics  of  the  flow  solver  NSMB  and  the  numerical  approach 
used  to  compute  time-accurate  solutions  around  a  moving  grid  are  first  presented.  A  few  words  are 
then  mentioned  about  the  geometry  chosen  and  the  mesh  generation  process.  The  main  part  of  the 
paper  deals  with  presentation  of  the  testcases  and  analysis  of  the  results  obtained.  This  part  is 
divided  between  inviscid,  laminar  and  turbulent  results.  Some  pictures  of  the  flowfield  are  then 
included  in  a  part  on  visualization.  Finally,  the  method  used  to  reduce  the  amount  of  data  stored  on 
disk  that  implies  the  co-processing  of  the  simulation  and  the  visualization  is  discussed  and 
illustrated  with  a  2D  example. 
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Numerical  Approach 

The  flow  solver  used  for  this  study  is  NSMB  (Navier-Stokes  Multi  Block)  developed  in  a  joint 
research  project  between  European  universities  and  aerospace  industries26.  NSMB  is  a  cell-centred 
finite  volume  CFD  code  that  solves  the  Navier-Stokes  equations  on  block-structured  grids.  Several 
turbulence  models  are  implemented  including  algebraic,  one-equation  and  two-equations  models 
and  it  is  possible  to  simulate  chemically  reacting  flows.  The  inviscid  flux  vectors  at  cell  faces  are 
approximated  using  either  Jameson's  central  scheme  or  a  second  order  Total  Variation 
Diminishing  (TVD)  version  of  Roe's  upwind  scheme.  For  the  central  scheme,  second-order 
artificial  dissipation  terms  are  added  near  discontinuities  and  fourth-order  terms  are  used  to 
suppress  odd/even  oscillations.  The  integration  in  time  of  the  resulting  equations  is  done  either  by 
the  explicit  Runge-Kutta  method  or  the  implicit  Lower-Upper-Symmetric  Gauss  Seidel  (LUS-GS) 
method. 

All  the  results  presented  in  this  paper  are  obtained  from  inviscid,  laminar  (viscous  non 
turbulent)  and  turbulent  calculations  done  with  the  2nd  order  upwind  scheme.  All  the  calculations 
were  done  with  a  serial  version  of  NSMB  running  on  a  3-processor  Fujitsu  VX  machine  with  a 
peak  performance  of  6.6  GFlops. 

Dual  Timestepping 

To  solve  the  time-accurate  flowfield  of  the  present  study,  the  dual  time-stepping  method  is 
used.  This  technique  introduces  an  outer  time-stepping  loop  for  a  real  time-accurate  timestep  using 
a  fully  implicit  scheme,  and  an  inner  loop  with  a  fictitious  timestep  to  reach  a  “steady  state”  at 
each  real  timestep.  The  advantage  of  the  method  is  that  acceleration  techniques  like  local  time¬ 
stepping,  multigridding  and  implicit  schemes  can  be  used  in  the  inner  loop  to  iterate  the  solution  to 
the  “steady  state”.  More  information  on  the  method  can  be  found  in  a  previous  paper17  by  the  same 
authors. 

The  dual  timestepping  method  can  be  presented  in  a  more  visual  way  when  one  looks  at  the 
convergence  history.  Figure  1  shows  the  L2  residuals  for  the  density  p  as  function  of  the  number 
of  iterations.  This  is  a  typical  choice  of  parameters  used  in  the  present  study:  60  outer  timesteps 
per  oscillation  and  100  inner  iterations.  For  each  outer  or  real  timestep,  one  clearly  sees  the 
residuals  decrease  during  the  100  inner  iterations  aiming  at  the  “steady”  state.  In  addition  to  the 
residuals,  the  evolution  of  the  angle  of  attack  is  plotted  with  the  scale  on  the  right.  It  is  interesting 
to  see  that  the  level  of  convergence  is  low  for  high  angles  of  attack  and  starts  to  improve  on  the 
downstroke  motion  when  a  becomes  less  than  10°  reaching  its  limit  on  the  upstroke  motion  around 
22°. 

Moving  Grid 

The  pitching  computations  are  performed  with  a  moving  grid  (the  so-called  ALE,  Arbitrary 
Lagrangian  Eulerian  approach).  The  whole  grid  rotates  like  a  rigid  body  with  the  pitched  object 
and  unsteady  boundary  conditions  are  implemented  for  the  far  field  and  surface  boundaries.  At  the 
wall,  the  velocity  is  equal  to  the  velocity  of  the  moving  wall  for  viscous  calculations  while  for 
inviscid  calculations,  only  the  normal  components  of  the  velocities  are  equal.  Reference17  contains 
some  more  information  with  the  governing  equations  that  take  into  account  the  mesh  velocity. 

Wing  Geometry  and  Mesh  Generation 

The  geometry  chosen  for  the  present  study  is  the  model  used  by  Soltani  and  Bragg20.  It  is  a  flat 
plate  70°-swept  wing  with  a  0.25  inch  (6.4  mm)  thickness  and  sharpened  edges  using  a  14°  bevel 
on  the  lowel  side,  see  Figure  2.  The  mesh  used  for  all  the  computations  presented  below  was  made 
by  S.  Gortz8,9  using  the  mesh  generator  ICEM  CFD  Hexa.  It  is  a  6-block  C-O  type  mesh  with  more 
than  700,000  volume  cells,  see  Figure  3.  The  chosen  topology  results  in  a  singularity  at  the  apex. 
To  resolve  the  flowfield  above  the  wing,  130  points  were  placed  in  the  chordwise  direction  along 
the  wing  and  40  points  in  the  spanwise  direction.  Only  half  of  the  wing  is  represented  in  the  mesh 
assuming  that  in  the  case  of  no  sideslip  the  longitudinal  middle  plane  is  a  symmetry  plane.  The  far- 
field  boundaries  extent  one  root  chord  upstream,  and  two  downstream  of  the  wing  and  one  root 
chord  in  the  normal  direction. 

(Several  other  meshes  are  currently  being  tried  on  the  same  geometry  but  without  any  success. 
They  are  H-H  type  grids.  The  first  one  was  used  by  Gortz8  in  his  first  studies  and  proved  to  give 
satisfactory  static  steady  results.  The  second  mesh  has  a  finer  structure  right  above  the  wing  to 
better  capture  the  vortical  flow.  The  computations  started  with  both  meshes  failed  to  converge  for 
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a  complete  oscillation  around  the  mean  angle  of  attack  22°.  That  is  the  reason  why  the  results  are 
not  included  in  the  present  paper.) 


Results  and  Analysis 

In  this  section,  the  results  of  the  computations  on  the  pitching  delta  wing  are  presented.  The 
pitching  motion  consists  of  sinusoidal  oscillations  defined  by  the  angle  of  attack  law: 
a(t)=am  +  Aasin{cot) 


Several  combinations  of  the  mean  angle  and  amplitude  have  been  computed  but  most  of  the 
effort  has  been  put  on  oscillations  around  a  =  22°  with  semi-amplitude  Aa  =  18°.  For  all  the 
computed  cases,  the  axis  of  rotation  is  placed  at  x/c  =  0.4  and  the  freestream  Mach  number  is 
Moo  =  0.2.  The  reduced  frequency  of  the  oscillations  is  the  same  for  all  the  computations: 
k  =0.0376  corresponding  to  one  of  the  test  cases  experimentally  studied  by  Brandon3.  This 
reduced  frequency  corresponds  to  an  oscillation  frequency  of  CO  =  5.12  rad/s  and  the  period  of  the 
oscillation  is  T  -  1.23  s.  For  each  calculation,  the  result  of  a  static  computation  at  the  mean  angle 
of  attack  is  used  as  starting  point  for  the  pitching  computation. 

Static  computations 

As  static  computations  are  required  before  starting  pitching  calculations,  few  static  steady 
results  are  available  even  if  it  is  not  the  main  point  of  interest  of  the  study  (a  more  complete 
analysis  of  static  steady  results  can  be  found  in  references  8  and  9  where  Gortz  presents  both 
Navier-Stokes  and  Euler  computations  for  the  same  geometry  and  using  the  same  flow  solver). 
Figure  4  shows  results  for  the  static  normal  force  coefficient  CN  compared  to  experimental  data3. 
In  this  picture,  all  the  numerical  results  are  obtained  with  the  second  order  upwind  scheme.  Good 
agreement  is  achieved  for  angles  between  10°  and  35°.  It  is  important  to  notice  that  the  turbulent 
results  (only  obtained  at  12°,  22°  and  38°)  are  very  close  to  the  inviscid  results.  Poor  agreement  is 
obtained  for  the  Euler  results  at  low  angles  of  attack.  However,  a  viscous  non-turbulent  calculation 
(“Laminar”  in  Figure  4)  shows  a  better  prediction  for  a  =  0°  indicating  that  the  origin  of  the 
problem  is  due  to  viscosity  effects  (a  turbulent  calculation  is  needed  to  confirm  this).  A  possible 
explanation  could  be  the  decambering  effect  of  the  boundary  layer  which  is  due  to  the  different 
thickness  of  the  boundary  layers  on  the  upper  and  lower  sides  of  the  wing,  that  slightly  changes  the 
camber  of  the  section  seen  by  the  flow. 

Euler  computations 

Main  testcase:  a{t)-22°  + 18 °  sin  (ox) 


The  experimental  testcase  chosen  to  compare  with  the  numerical  results  is  the  wind  tunnel 
experiment  by  Brandon3:  pitching  oscillations  of  a  70°  delta  wing  around  an  angle  of  attack  of  22° 
and  with  a  semi-amplitude  of  18°.  Unfortunately  only  the  normal  force  coefficient  CN  is  reported 
in  Figure  3  of  reference  3  and  no  data  is  available  for  the  drag  or  pitching  moment  coefficients. 

The  results  for  the  normal  force  coefficient  in  Figure  5  reproduce  the  expected  hysteresis  loops 
as  the  values  are  higher  for  the  upstroke  motion  than  for  the  downstroke  as  observed  in  the 
experiment.  The  maximal  normal  force  is  higher  than  the  experimental  values  and  this  is  true  for 
all  the  computations  performed  in  the  present  study.  A  possible  cause  of  this  difference  could  be 
the  experimental  mounting  system  that  has  a  major  influence  on  the  flow  around  the  wing  (see  for 
example  Ericsson4).  However,  the  discrepancy  observed  for  low  angles  of  attack  (a<15°)  is 
corrected  by  taking  into  account  the  fluid  viscosity  as  explained  earlier  for  the  static  results  (see 
also  the  part  on  the  laminar  calculations). 

(The  results  for  the  drag  and  pitching  moment  coefficients  are  included  in  Figure  6,  see  next 
paragraph). 

Other  testcases 

Other  combinations  of  the  mean  angle  of  attack  and  semi-amplitude  have  been  computed  so 
that  results  for  the  following  motions  are  available. 


•  Case  1:  large  am,  large  Aa: 

•  Case  2:  large  am,  low  Aa: 

•  Case  3:  low  am,  large  Aa: 

•  Case  4:  low  am,  low  Aa: 


ct(t>)—22  +  18  sin(cot) 
a(t)-38°  +6°  sin  (ox) 
cc{t^— 12  +18  sin(cot) 
a(t)-12°  +3°  sin(cot) 
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The  results  are  reported  in  Figure  6,  inspired  of  Fig. 6  in  Hummel  and  Loeser11.  The  results  for 
case  4  are  always  aligned  and  no  hysteresis  can  be  seen  because  the  vortex  breakdown  is  not 
located  above  the  wing  for  this  range  of  angles  of  attack.  For  the  normal  force  coefficient  in  Figure 
6(a),  the  results  for  case  3  extend  the  hysteresis  from  case  1  towards  negative  angles  of  attack 
following  the  same  trend.  It  is  to  mention  that  the  hysteresis  loop  for  case  3  is  quite  thin  and 
bended  around  10°.  In  fact,  as  the  angle  of  attack  increases,  the  width  of  the  hysteresis  loops 
increases  too  and  the  loops  slightly  rotate  in  the  clockwise  direction.  This  is  even  more  noticeable 
in  Hummel  and  Loeser's  figure  where  the  loops  around  42°  and  48°  are  horizontal. 

The  drag  results  in  Figure  6(b)  are  also  in  good  qualitative  agreement  with  Hummel  and 
Loeser's  results:  a  higher  drag  is  obtained  in  the  upstroke  motion  creating  hysteresis  loops.  The 
dynamic  results  follow  the  static  curve  and  only  deviate  slightly  at  high  angles  of  attack.  The 
pitching  moment  coefficient  in  Figure  6(c)  shows  almost  a  linear  behaviour  with  very  thin  and  flat 
hysteresis  loops  only  diverging  from  the  linear  trend  at  large  angles  of  attack.  This  figure  is 
different  from  Hummel  and  Loeser's.  Their  results  show  wide  and  large  hysteresis  loops  around 
the  static  curve  whereas  in  our  case  the  difference  between  upstroke  and  downstroke  pitching 
moment  is  small.  This  is  certainly  due  to  the  higher  reduced  frequency  used  in  the  experiment 
( k  =0.28)  which  is  known  to  create  larger  over-  and  -undershoot  in  the  dynamic  results  compared 
to  static  values. 

Laminar  computations 

The  major  part  of  this  study  is  concerned  with  in  viscid  calculations  but  a  few  laminar 
computations  (viscous  non-turbulent)  have  also  been  performed  for  a  Reynolds  number  Re  equal 
to  3.16E6.  As  already  seen  for  the  static  case,  the  viscosity  has  an  influence  on  the  results  at  low  to 
moderate  angles  of  attack.  Figure  7  shows  that  the  laminar  normal  force  coefficient  CN  is  closer  to 
the  experimental  values  for  a<l  5°,  considerably  reducing  the  gap  observed  between  in  viscid  and 
experimental  results.  The  possible  cause  of  this  better  agreement  (decambering  effect  of  the 
boundary  layer)  has  been  presented  in  the  section  on  the  static  results.  On  the  other  hand,  for 
angles  of  attack  higher  than  20°  no  difference  is  noticed  between  inviscid  and  laminar  results.  For 
the  drag  coefficient  results,  no  difference  between  laminar  and  inviscid  calculations  is  noticeable 
either. 

Turbulent  computations 

Once  the  inviscid  calculations  prove  to  give  reasonable  results,  the  next  step  is  to  try  turbulent 
calculations  and  see  if  the  addition  of  viscosity  and  turbulence  improves  the  accuracy  or  not.  The 
main  testcase  (case  1)  has  been  tried  first  with  both  the  algebraic  model  of  Baldwin  and  Lomax 
and  the  one-equation  model  of  Spalart  and  Allmaras.  The  results  for  the  normal  force  coefficient 
are  plotted  in  Figure  8  with  the  inviscid  results  and  the  experimental  ones.  Unfortunately,  at  the 
time  of  the  writing  of  this  paper,  the  computations  have  not  successfully  converged  for  a  complete 
oscillation  which  explains  why  only  half  of  the  cycle  is  plotted  for  the  turbulent  data.  The  first 
important  remark  is  that  for  the  high  angles  of  attack  the  turbulent  results  only  slightly  differ  from 
the  inviscid  ones:  the  maximum  normal  force  is  a  bit  lower  in  the  turbulent  case  and  the  same 
applies  to  the  first  half  of  the  oscillation  in  general  but  the  difference  is  very  small.  No  difference 
can  be  noticed  between  the  algebraic  and  the  one-equation  model  except  maybe  for  the  lower 
angles  of  attack  where  the  one-equation  model  seems  to  follow  the  inviscid  curve  whereas  the 
algebraic  model  yields  slightly  higher  values.  Of  course  a  complete  oscillation  would  be 
interesting  to  confirm  the  trend  and  see  if  the  turbulent  results  behave  like  the  laminar  ones  or  not. 

Case  2  with  oscillations  around  38°,  which  is  a  less  demanding  test  case  than  case  1  as  the 
amplitude  of  the  oscillation  is  rather  small,  has  been  successfully  computed  with  the  one-equation 
turbulent  model  of  Spalart  and  Allmaras.  The  normal  force  results  are  presented  in  Figure  9  and 
compared  to  the  inviscid  values.  The  general  shape  of  the  hysteresis  loop  is  the  same  for  both 
computations  although  the  turbulent  results  are  slightly  higher  in  the  downstroke  part  of  the 
oscillation.  One  can  also  notice  that  the  turbulent  loop  is  smoother  and  rounder  (see  the  upstroke 
part  where  the  inviscid  results  coincide  with  the  turbulent  values  except  for  some  points  out  of  the 
oval  shape  of  the  loop).  Finally,  the  static  results  at  the  mean  angle  of  attack  of  38°,  also  plotted  in 
Figure  9,  are  almost  similar.  The  results  for  the  drag  and  the  pitching  moment  coefficient  (not 
reported  in  the  paper)  show  exactly  the  same  behaviour  and  all  the  remarks  made  for  the  normal 
force  coefficient  apply  to  the  other  two  coefficients  too. 
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The  work  on  the  turbulent  simulations  is  going  on  and  hopefully  more  results  are  soon 
available. 

Flowfield  visualization 

The  visualization  of  the  unsteady  flowfield  around  the  pitching  delta  wing  has  not  been  so  far 
the  major  point  of  the  study  and  only  few  pictures  are  available.  However,  Figure  10  shows  some 
views  of  the  vortex  above  the  wing  at  different  angles  of  attack  for  the  main  testcase  (case  1).  The 
wing  is  coloured  by  pressure  and  the  zone  of  low  pressure  is  visible  below  the  vortex.  The  vortex 
is  made  visible  by  creating  an  isosurface  of  total  pressure. 

Figure  10(a)  shows  the  vortex  for  a  static  steady  calculation  at  the  mean  angle  of  attack  22° 
while  Figure  10(b)  to  Figure  10(i)  show  the  behaviour  of  the  vortex  for  a  complete  oscillation  from 
an  angle  of  attack  of  4°,  going  up  to  40°  and  then  down  again  to  4°.  The  formation  and 
disappearance  of  the  vortex  is  well  reproduced  as  well  as  the  motion  of  the  breakdown  that  can  be 
associated  with  the  sudden  increase  of  the  section  of  the  isosurface  (particularly  visible  on  Figure 
10(e)  and  Figure  10(h)  at  the  front  part  of  the  wing).  It  is  interesting  to  qualitatively  compare  the 
shape  of  the  vortex  and  the  position  of  the  breakdown  at  a  constant  angle  of  attack.  Figure  10(e) 
and  Figure  10(g)  are  both  taken  at  an  angle  of  attack  of  34.73°,  but  the  shape  of  the  vortex  is 
different:  in  the  upstroke  motion  (Figure  10(e)),  the  vortex  is  not  burst  at  the  apex  of  the  wing  but 
slightly  downstream,  and  the  expansion  of  the  core  is  clear;  in  the  downstroke  motion  (Figure 
10(g))  on  the  contrary,  the  core  is  already  expanded  at  the  apex  of  the  wing.  This  is  expected  as  the 
breakdown  in  the  latter  case  is  moving  downstream  from  its  position  at  the  apex  at  40°.  The  same 
comparison  can  be  made  for  the  pictures  at  the  mean  angle  of  attack  22°  although  the  difference 
between  the  static  (Figure  10(a))  and  the  upstroke  (Figure  10(d))  cases  is  not  so  clear.  At  4°,  a 
small  zone  of  lower  pressure  is  seen  on  the  wing  but  the  vortex  structure  is  to  weak  and  not  seen 
for  this  level  of  total  pressure.  All  the  pictures  in  Figure  10  are  generated  from  an  output  file  that 
contains  13  solutions  for  the  whole  flowfield  at  13  different  angles  of  attack.  The  size  of  such  a  file 
is  rather  large  (around  850  Mb)  and  deters  from  saving  the  solution  at  each  time  step,  although  it 
would  be  necessary  to  see  the  formation  and  burst  of  the  vortex  in  a  more  accurate  way.  That  is 
why  the  new  saving  method  presented  below  is  introduced. 

Co-processing  simulation  and  visualization  using  flow  feature  extraction 

The  work  presented  in  this  part  is  inspired  by  a  paper  by  Haimes  and  Jordan10  who  present  a 
new  approach  to  analyse  3D  large  scale  transient  solutions  from  CFD  simulations.  In  the  paper,  the 
classical  way  of  processing  CFD  results  i.e.  solving  the  problem  with  the  flow  solver  in  a  first  step 
and  then  analyse/visualize  the  results  in  a  post-processing  step,  is  shown  to  be  unsuitable  for  large 
testcases  due  to  the  amount  of  data  necessary  to  store  on  disk.  This  is  particularly  true  for  transient 
simulations  where  the  solution  need  to  be  saved  at  several  time  steps.  Haimes  and  Jordan  rather 
suggest  to  save  only  the  “interesting”  part  of  the  solution  to  the  disk.  This  implies  that  the  analysis 
of  the  solution  is  done  prior  to  any  storage  on  disk.  In  order  to  do  that,  the  flow  solver  and  the 
visualization  tool  need  to  be  running  and  communicating  simultaneously:  as  soon  as  a  solution  is 
obtained  the  visualization  tool  analyses  it  and  extracts  the  useful  information  that  is  the  only  part 
saved  to  disk.  To  fully  take  advantage  of  the  method,  the  analysis  of  the  results  has  to  be  made 
automatic  so  as  to  avoid  that  the  flow  solver  waits  for  the  user  to  select  what  to  save  to  disk.  This 
so-called  “interactive  visualization”  has  the  drawback  of  necessitating,  in  advance,  quite  a  good 
knowledge  of  the  future  solution.  For  example,  the  user  must  know  beforehand  if  the  solution  may 
contain  vortices  in  order  to  have  the  visualization  tool  search  for  them  in  the  flowfield.  This  is 
compared  by  Haimes  and  Jordan  to  experimental  tests  where  a  good  idea  of  the  regions  of  interest 
is  required  before  positioning  the  probes.  Examples  of  flow  features  that  can  be  extracted  are 
vortex  cores  or  shocks  that  requires  special  extracting  algorithms  but  also  cutting  planes  or  simply 
parts  of  the  flowfield  like  walls  can  be  thought  of.  For  these  parts  only,  the  same  information  as  in 
the  classical  approach  is  stored  for  future  probing  and  analysis,  which  results  in  a  several  orders  of 
magnitude  lower  storage  requirement.  Of  course  it  is  still  possible  to  save  the  entire  flowfield  to 
disk  but  more  seldom  i.e.  every  n  time  steps  instead  of  each  time  step  for  an  unsteady  simulation. 
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The  interactive  visualization  is  of  course  very  attractive  when  computing  a  pitching  delta  wing 
which  is  both  a  large  scale  3D  testcase  and  a  transient  one.  That  is  why  this  method  has  been 
investigated  to  help  reducing  the  amount  of  data  generated  in  the  present  study  (and  future  large 
unsteady  simulations).  The  basics  of  the  method  are  illustrated  by  the  diagram  in  Figure  11:  the 
communication  between  the  flow  solver  NSMB  and  the  visualization  software  EnSight  is  done 
through  the  MEMCOM  database  which  is  the  output  file  from  NSMB.  This  is  obviously  not  the 
most  efficient  communication  means  due  to  writing  and  reading  time  losses  but  it  has  the 
advantage  of  being  quite  simple  to  implement  when  both  the  flow  solver  and  the  visualization  tool 
are  given  as  finished  products,  not  necessarily  designed  to  communicate  with  each  other.  In  the 
present  case,  NSMB  writes  a  solution  to  the  MEMCOM  database  (in  the  usual  way)  and  an 
additional  routine  calls  EnSight  to  process  the  data.  The  input  to  EnSight  is  a  command  file 
usually  created  while  visualizing  the  static  solution  that  is  necessary  as  starting  point  to  the 
unsteady  moving  grid  calculation.  At  the  end  of  the  data  processing,  the  visualization  software 
saves  to  disk  the  interesting  geometric  entities  with  the  variables  the  user  wants  to  save.  NSMB 
then  continues  with  the  computation  at  the  next  time  step.  The  solution  saved  in  the  database  is 
either  overwritten  by  the  solution  at  the  next  time  step  or  kept  in  the  file  for  a  permanent  storage 
(and  future  post-processing  of  the  entire  flowfield). 

Example  of  a  pitching  NACA  0012  airfoil: 

At  the  time  of  the  writing  of  this  paper,  no  example  for  the  pitching  delta  wing  is  available  yet 
but  a  2D  validating  testcase  has  been  successfully  run:  an  Euler  calculation  of  the  NACA0012 
airfoil  pitching  around  a  mean  angle  of  attack  of  5°  and  with  a  semi  amplitude  of  5°.  As  NSMB  is  a 
3D  solver,  the  mesh  consists  of  very  closely  spaced  planes  orthogonal  to  the  airfoil  (4  in  total  that 
build  a  span  of  3  cells)  and  in  all  there  are  8192  cells  in  the  C-type  mesh.  The  computation 
methods  used  are  the  same  as  in  the  delta  wing  case:  ALE  and  dual  time  stepping  with  40  inner 
steps  and  40  outer  in  one  oscillation.  The  MEMCOM  database  containing  1 1  solutions  (that  is  1 1 
sets  of:  density,  temperature,  Mach  number,  pressure,  velocity  field  and  residuals)  saved  every  4th 
outer  time  steps  and  including  the  mesh  information  has  a  size  of  13Mb.  Saving  only  the  geometry 
for  the  airfoil  surface  and  a  longitudinal  plane  of  the  grid  with  the  pressure  and  Mach  number 
information  reduces  the  disk  storage  to  960Kb  (and  even  440Kb  if  only  one  geometry  is  saved  as 
is  the  case  with  the  MEMCOM  database)  which  is  only  7.3%  of  the  original  data  size.  It  is  then 
possible  to  produce  pictures  like  Figure  12  (visualization  of  the  Mach  number  in  the  longitudinal 
plane  and  pressure  on  the  airfoil)  every  4th  time  step.  For  a  smaller  disk  storage  than  in  the  original 
case  with  the  MEMCOM  database,  it  could  even  be  possible  to  save  the  solution  at  each  timestep. 

Conclusion  and  outlook 

In  the  present  paper,  the  results  from  Euler  and  Navier-Stokes  simulations  of  a  70°-swept  delta 
wing  have  been  presented.  The  main  effort  is  put  on  predicting  the  hysteresis  loops  that  form  in  the 
aerodynamic  coefficients  history.  Both  the  inviscid  and  turbulent  calculations  successfully 
reproduce  the  hysteresis  effect  without  any  major  difference  between  them  except  at  low  angles  of 
attack.  However,  the  lack  of  turbulent  results  for  these  angles  prevent  us  from  drawing  any 
conclusion.  The  prediction  of  the  maximal  normal  force  coefficient  is  not  improved  by  turbulence: 
the  maximal  normal  force  is  almost  the  same  for  turbulent  and  inviscid  calculations  and  in  both 
cases,  higher  than  the  experimental  value. 

The  flowfield  visualization  shows  the  formation  and  disappearance  of  the  vortex  above  the 
wing.  Even  if  a  picture  of  the  flowfield  is  only  available  at  some  points  in  the  oscillation,  the 
results  are  very  promising  and  encouraging  for  a  more  complete  and  time-accurate  exploration  of 
the  flowfield. 

Finally,  the  introduction  of  the  co-processing  of  the  flow  solver  and  the  visualization  tool 
proves  to  highly  reduce  the  amount  of  data  stored  on  disk.  If  the  application  of  this  method  to  the 
case  of  the  pitching  delta  wing  is  successful,  the  future  computations  and  visualizations  should  be 
more  accurate  in  time  as  the  timesteps  could  be  reduced  (and  more  solutions  saved).  This  is  the 
next  step  in  this  study. 

Another  important  matter  is  of  course  to  have  the  turbulent  computations  converge  and  get  a 
complete  oscillation  for  the  main  testcase  with  oscillations  around  22°.  Some  other  future  tasks  are 
to  try  a  different  mesh  around  the  same  geometry  and  see  the  influence  of  the  fineness  or  topology 
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of  the  grid  on  the  unsteady  results.  Finally,  the  author  also  intends  to  try  a  4th  order  accurate  central 
scheme  that  yields  promising  static  results. 
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Convergence-  and  a-history  Static  results 


Figure  1:  Convergence  and  a  history  for  a 
calculation  with  60  outersteps  per 
oscillation  and  100  inner  iterations 
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Figure  2:  Wing  geometry 


C-0  mesh 


Figure  4:  Static  normal  force  results 
compared  to  experimental  data3 


Normal  force  coefficient 


Figure  5:  Normal  force  results  from  inviscid 
calculations  compared  to  experiment3 


Static  and  dynamic  normal  force  coefficient 


(a)  Normal  force  coefficient 


Figure  3:  C-O  mesh 
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Static  and  dynamic  drag  coefficient 


(b)  Drag  coefficient 


Static  and  dynamic  pitching  moment  coefficient 


(c)  Pitching  moment  around  the  rotation 
axis  at  x/c  =  0.4 

Figure  6:  Static  and  dynamic  results  for 
various  mean  angles  of  attack  and 
amplitudes 


Influence  of  viscosity  on  normal  force  coefficient 


Preliminary  turbulent  results 


Figure  8:  Preliminary  normal  force 
coefficient  results  from  turbulent 
computations  for  case  1 

Inviscid  and  turbulent  (Spalart-Allmaras  model)  results,  oscillations  38°+/-6° 


Figure  9:  Normal  force  coefficient  results  for 
a  turbulent  calculation  (case  2) 


(a)  Static  a  =  22° 


Figure  7:  Normal  force  coefficient  results  for 
a  laminar  computation 
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Figure  10:  Flow  visualization:  vortex  visualized  by  an  isosurface  of  total  pressure  (case  1) 
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Figure  1 1 :  Interactive  visualization 
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Figure  12:  NACA-0012  visualization 
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